function G_dB = load_vehicular_option2(fr, location)
%LOAD_VEHICULAR_OPTION2(fr, location)
% Loads option 2 antenna power patterns for vehicular UEs from 
% 3GPP TR 37.885 V15.1.0 (2018-09)
%
% Input:
%   fr          frequency range
%               1: below 6GHz
%               2: above 6GHz
%   location    location of antenna on vehicle
%               1: front bumper
%               2: front rooftop
%               3: rear rooftop
%               4: rear bumper
% 
% Output:
%   G_dB        antenna power pattern
%
%
% QuaDRiGa Copyright (C) 2011-2019
% Fraunhofer-Gesellschaft zur Foerderung der angewandten Forschung e.V. acting on behalf of its
% Fraunhofer Heinrich Hertz Institute, Einsteinufer 37, 10587 Berlin, Germany
% All rights reserved.
%
% e-mail: quadriga@hhi.fraunhofer.de
%
% This file is part of QuaDRiGa.
%
% The Quadriga software is provided by Fraunhofer on behalf of the copyright holders and
% contributors "AS IS" and WITHOUT ANY EXPRESS OR IMPLIED WARRANTIES, including but not limited to
% the implied warranties of merchantability and fitness for a particular purpose.
%
% You can redistribute it and/or modify QuaDRiGa under the terms of the Software License for 
% The QuaDRiGa Channel Model. You should have received a copy of the Software License for The
% QuaDRiGa Channel Model along with QuaDRiGa. If not, see <http://quadriga-channel-model.de/>. 


if ~( any( numel(fr) == 1 ) && isnumeric(fr) && isreal(fr) &&...
        max(fr)<=2 && min(fr)>=1 )
    error('QuaDRiGa:qd_arrayant:wrongInputValue','??? "fr" must be a numbel between 1 and 2')
end

if ~( any( numel(location) == 1 ) && isnumeric(location) && isreal(location) &&...
        max(location)<=4 && min(location)>=1 )
    error('QuaDRiGa:qd_arrayant:wrongInputValue','??? "location" must be a number between 1 and 4')
end

switch fr
    case 1
        switch location
            case 1
                %% front bumper from Table 6.1.4-10A: Front bumper antenna element pattern for vehicle UE in Option 2 for 6 GHz
                abc_G_phi1 = [...
                    27.855	0.9281	-1.3062; ...
                    31.967	2.7336	1.4674; ...
                    23.665	3.0514	-1.8167; ...
                    1.8441	7.7354	-0.6126; ...
                    ];
                abc_G_phi2 = [...
                    146.47	1.095	-1.6043; ...
                    131.87	1.192	1.4566; ...
                    5.2383	3.9996	-0.8389; ...
                    1.9429	9.0248	1.6705; ...
                    ];
                
                ab_G_theta2 = [...
                    -20.793	0; ...
                    -3.4026	0.0721; ...
                    1.3366	0.1213; ...
                    2.1764	-0.0595; ...
                    0.5881	0.0102; ...
                    -2.4068	0.0675; ...
                    1.7561	-0.1513; ...
                    0.7899	-0.0133; ...
                    -1.0144	0.0177; ...
                    ];
                
                gain_dBi = 13;
                
            case 2
                %% front rooftop from Table 6.1.4-10B: Front rooftop antenna element pattern for vehicle UE in Option 2 for 6 GHz
                abc_G_phi1 = [...
                    27.156	0.8354	-0.9866; ...
                    10.637	1.2651	1.9323; ...
                    2.1018	3.7081	-1.6452; ...
                    0.9267	12.057	-2.9006; ...
                    0.8573	8.5269	3.5787; ...
                    0.7456	15.27	-2.9591; ...
                    ];
                abc_G_phi2 = [...
                    18.648	1.0206	-0.9611; ...
                    11.099	2.2755	1.7808; ...
                    213.68	3.7367	-1.6448; ...
                    208.11	3.7558	1.492; ...
                    ];
                
                ab_G_theta2 = [...
                    -11.719	0; ...
                    -2.4504	0.0053; ...
                    -0.9825	-0.0508; ...
                    -0.0419	-0.0225; ...
                    0.2384	0.0065; ...
                    0.2406	-0.0023; ...
                    0.1532	0.006; ...
                    0.1285	0.0062; ...
                    0.1407	0.0038; ...
                    ];
                
                gain_dBi = 3;
            case 3
                %% rear rooftop from Table 6.1.4-10C: Rear rooftop antenna element pattern for vehicle UE in Option 2 for 6 GHz
                abc_G_phi1 = [...
                    21.785	0.7541	-2.5975; ...
                    122	2.0393	-0.3421; ...
                    4.5185	4.4006	2.4657; ...
                    25.257	3.2885	-0.6558; ...
                    125.16	2.3309	2.7194; ...
                    ];
                abc_G_phi2 = [...
                    89.841	0.5245	-0.2707; ...
                    64.942	0.6058	2.9391; ...
                    2.4779	4.1726	-1.4849; ...
                    1.097	7.9882	-1.7998; ...
                    ];
                ab_G_theta2 = [...
                    -12.399	0; ...
                    -0.6788	-0.0034; ...
                    -1.025	-0.0364; ...
                    -0.0755	-0.0404; ...
                    0.4216	0.014; ...
                    -0.6139	-0.0053; ...
                    0.2297	0.0063; ...
                    -0.4418	-0.0038; ...
                    0.0055	-0.0141; ...
                    ];
                gain_dBi = 3;
                
            case 4
                %% rear bumper from Table 6.1.4-10D: Rear bumper antenna element pattern for vehicle UE in Option 2 for 6 GHz
                abc_G_phi1 = [...
                    28.64	0.6822	-1.6819; ...
                    17.941	1.0841	0.8466; ...
                    10.034	4.5559	-2.317; ...
                    9.9931	4.9505	0.8248; ...
                    103.05	9.2847	0.5153; ...
                    102.49	9.3113	-2.6065; ...
                    ];
                abc_G_phi2 = [...
                    29.045	0.2614	-0.7172; ...
                    7.9836	2.8649	2.965; ...
                    5.1881	3.3343	-0.1509; ...
                    0.4384	7.4621	1.8565; ...
                    1.0833	47.882	2.0727; ...
                    1.131	35.828	-0.901; ...
                    1.0114	33.576	-0.1873; ...
                    1.0014	39.562	-1.0305; ...
                    ];
                ab_G_theta2 = [...
                    -17.964	0; ...
                    4.1952	0.0073; ...
                    -0.3746	-0.0202; ...
                    1.2325	0.005; ...
                    0.2862	-0.0764; ...
                    -0.2796	0.0252; ...
                    -0.063	-0.0114; ...
                    -1.0205	0.0298; ...
                    -0.2964	0.0095; ...
                    ];
                gain_dBi = 11;
        end
    case 2
        switch location
            case 1
                % 30 and 63 GHz
                %% front bumper from Table 6.1.4-11A: Front bumper antenna element pattern for vehicle UE in Option 2 for 30 and 63 GHz
                abc_G_phi1 = [...
                    103.4	1.2756	-1.4588; ...
                    93.054	1.4402	1.6919; ...
                    1.2433	16.391	-0.2701; ...
                    ];
                
                abc_G_phi2 = [...
                    97.485	1.2777	-1.5428; ...
                    88.914	1.4458	1.5626; ...
                    1.6315	6.6824	-1.4228; ...
                    ];
                
                ab_G_theta2 = [...
                    -14.712	0; ...
                    -0.814	-0.0126; ...
                    4.5025	-0.4615; ...
                    0.6572	-0.0404; ...
                    0.3252	0.5985; ...
                    -2.0187	0.0372; ...
                    -1.6555	0.2589; ...
                    -0.843	0.6858; ...
                    -0.2954	-0.2211; ...
                    ];
                
                gain_dBi = 16.6;
            case 2
                %% front rooftop from Table 6.1.4-11B: Front rooftop antenna element pattern for vehicle UE in Option 2 for 30 and 63 GHz
                abc_G_phi1 = [...
                    53.19	1.0466	-0.9824; ...
                    171.45	2.3513	2.2593; ...
                    154.83	2.4609	-0.8176; ...
                    2.166	6.2024	0.5122; ...
                    1.2277	9.599	0.5604; ...
                    1.1363	14.311	-0.4448; ...
                    ];
                
                
                abc_G_phi2 = [...
                    34.284	1.0869	-1.1043; ...
                    36.721	2.2823	1.4488; ...
                    4.2902	5.2423	1.1692; ...
                    1.7681	7.7706	-3.641; ...
                    14.094	2.7398	-1.5731; ...
                    1.2267	10.465	2.0934; ...
                    1.0796	14.285	2.8759; ...
                    0.8846	31.578	-1.1871; ...
                    ];
                ab_G_theta2 = [...
                    -16.085	0; ...
                    -7.903	0.1101; ...
                    -0.2972	-1.4639; ...
                    0.4083	-0.0866; ...
                    0.6706	1.0598; ...
                    -0.9804	-2.5714; ...
                    -0.9085	0.1326; ...
                    0.1433	-0.1372; ...
                    -0.1282	-0.2697; ...
                    ];
                gain_dBi = 11.5;
            case 3
                %% rear rooftop from Table 6.1.4-11C: Rear rooftop antenna element pattern for vehicle UE in Option 2 for 30 and 63 GHz
                abc_G_phi1 = [...
                    55.839	1.0927	-1.1273; ...
                    95.402	2.0245	2.1133; ...
                    2.9533	6.0111	1.5482; ...
                    1.5696	12.462	3.2683; ...
                    1.2073	10.016	2.2276; ...
                    70.244	2.1361	-0.8181; ...
                    1.4093	13.127	0.6659; ...
                    ];
                
                
                abc_G_phi2 = [...
                    42.818	0.9167	-0.9757; ...
                    18.065	1.895	1.5513; ...
                    8.8444	3.6653	-0.1862; ...
                    9.5301	7.6263	0.1674; ...
                    15.018	9.0456	0.2995; ...
                    20.098	8.565	-2.8308; ...
                    1.5077	12.14	-1.3725; ...
                    ];
                
                ab_G_theta2 = [...
                    -18.471	0; ...
                    3.6106	-0.1824; ...
                    -1.4888	1.3614; ...
                    0.5772	0.5331; ...
                    1.7125	1.3205; ...
                    -0.0573	-2.5904; ...
                    0.0628	-1.0024; ...
                    -1.0526	-1.2366; ...
                    0.0433	-0.1443; ...
                    ];
                
                gain_dBi = 12.5;
            case 4
                %% rear bumper from Table 6.1.4-11D: Rear bumper antenna element pattern for vehicle UE in Option 2 for 30 and 63 GHz
                abc_G_phi1 = [...
                    87.458	0.0103	-2.9955; ...
                    4.1908	2.1185	1.8618; ...
                    1.0935	127.65	2.3653; ...
                    1.2773	16.837	1.1214; ...
                    ];
                
                
                abc_G_phi2 = [...
                    17.631	0.823	-1.5005; ...
                    12.155	1.9599	1.2643; ...
                    15.226	5.6698	2.2669; ...
                    15.495	5.7995	-0.9411; ...
                    0.3371	10.294	-3.8515; ...
                    0.5965	11.41	0.4383; ...
                    ];
                
                ab_G_theta2 = [...
                    -12.742	0; ...
                    -1.2007	0.2451; ...
                    4.1916	-0.5651; ...
                    0.1901	0.453; ...
                    0.4235	0.6884; ...
                    -1.3848	-0.4145; ...
                    -0.2381	0.2875; ...
                    0.4488	0.5441; ...
                    0.1606	-0.2555; ...
                    ];
                gain_dBi = 13.6;
        end
end

%%

phi = (0:360)/180*pi;
theta = (0:180)/180*pi;

G_phi1 = sum((abc_G_phi1(:, 1)*ones(size(theta)).*sin(abc_G_phi1(:, 2)*(pi/2-theta)+abc_G_phi1(:, 3)*ones(size(theta)))), 1);
G_phi2 = sum((abc_G_phi2(:, 1)*ones(size(theta)).*sin(abc_G_phi2(:, 2)*(pi/2-theta)+abc_G_phi2(:, 3)*ones(size(theta)))), 1);
k = 0:numel(ab_G_theta2(:, 1))-1;
G_theta2 = sum((ab_G_theta2(:, 1)*ones(size(phi))).*cos(k'*phi) + (ab_G_theta2(:, 2)*ones(size(phi))).*sin(k'*phi), 1);

o2n = 1:90;
o2n2 = o2n.*(90-o2n);
o2h = 1:180;
o2h2 = o2h.*(180-o2h);

G_theta1_up = (G_phi1(1) + G_phi2(1))/2;
G_theta1_low = (G_phi1(end) + G_phi2(end))/2;

G_dB = (4*([0, o2n2, o2n2].*G_phi2).'*[0, o2h, 180 - o2h] + ...
    4*([0, o2n2, o2n2].*G_phi1).'*[180, 180 - o2h, o2h] + ...
    [0, o2n, 90-o2n].'*(G_theta2.*[0, o2h2, o2h2]) + ...
    [G_theta1_up*[90, 90-o2n], G_theta1_low*o2n].'*[0, o2h2, o2h2] ...
    )./(180*4*[0, o2n2, o2n2].'*ones(1, 361) + ...
    ones(181, 1)*90*[0, o2h2, o2h2]);

G_dB(end, :) = G_theta1_low;
G_dB(1, :) = G_theta1_up;
G_dB(:, 1) = G_phi1;
G_dB(:, end) = G_phi1;
G_dB(:, 181) = G_phi2;
G_dB = G_dB - max(G_dB(:)) + gain_dBi;
G_dB = flipud(G_dB); % change from zenith angle to elevation angle

end
